function [backgroundimg, bprocess] = CreateBackgroundImageAuto(nstartframe, nendframe, nminstepsize, fvideo, fsubtract, handles)

flen = length(fvideo);
fout = fvideo(1:flen-4);
fout = strcat(fout, 'background.mat'); 

load(fsubtract);

cc = [cr1, cr2, cg1, cg2, cb1, cb2, ch1, ch2, cs1, cs2, cv1, cv2];

v = mmreader(fvideo);

img = read(v,2);
sz = size(img);

xlen = sz(1);
ylen = sz(2);

NumberOfFrames = get(v,'NumberOfFrames');

if nstartframe == -1
    nstartframe = 1;
end

if nendframe == -1
    nendframe = NumberOfFrames;
end


ntotalimages = 50; %how many images to average background across
ndeltatotalimages = 100; %if ntotalimages is not enough, increase ntotalimages by this value
nmaxtotalimages = 300;  %maximum number of images to average across

minstfill = 5;
nstfill = 15;

if nstfill < minstfill
    nstfill = minstfill;
end


stopprocess = 0;
bprocess = 1;

while bprocess == 1 & stopprocess == 0 & ntotalimages <= nmaxtotalimages
    
    stepsize = round((nendframe-nstartframe)/ntotalimages);
    
    if stepsize < nminstepsize
        stepsize = nminstepsize;
        stopprocess = 1;
    end
    
    [meanpixels, stdpixels, meanobjs, stdobjs] = GetMeanBWPixels(v, nstartframe, nendframe, stepsize, cc, nstremove, nstfill, nbwarea, dthresh, xlen, ylen);        
    
    nmin = round(meanpixels+stdpixels);

    [bprocess, imbuf, countbuf] = CreateBackground(v, nstartframe, nendframe, stepsize, cc, nstremove, nstfill, nbwarea, dthresh, xlen, ylen, nmin);
           
    ntotalimages = ntotalimages + ndeltatotalimages;

end

imbuf(:,:,1) = imbuf(:,:,1)./countbuf;
imbuf(:,:,2) = imbuf(:,:,2)./countbuf;
imbuf(:,:,3) = imbuf(:,:,3)./countbuf;   

backgroundimg = uint8(imbuf);





function [bprocess, imbuf, countbuf] = CreateBackground(v, startframe, endframe, stepsize, cc, nstremove, nstfill, nbwarea, dthresh, xlen, ylen, nmin)

i = startframe;
endframe = endframe-stepsize;

bwall = ones(xlen,ylen);
countbuf = zeros(xlen,ylen);
imbuf = zeros(xlen,ylen,3);

bprocess = 1;

while i <= endframe 
    i1 = i;
    i2 = i+stepsize;
    
    img1 = read(v,i1);
    img2 = read(v,i2);
    
    dimg1 = double(img1);    
    dimg2 = double(img2);    
    
    [bw1, idif1] = ProcessImage(img1, double(img2), double(rgb2hsv(img2)), cc, nstremove, nstfill, nbwarea, dthresh);
    [bw2, idif2] = ProcessImage(img2, double(img1),double(rgb2hsv(img1)), cc, nstremove, nstfill, nbwarea, dthresh);
     
    bw = bw1+bw2;   
    
    bwsum = sum(sum(bw));
     
    if bwsum > nmin
        for j = 1:xlen
            for k = 1:ylen
                if bw(j,k) == 0
                    bwall(j, k) = 0;
                    imbuf(j,k,1) = imbuf(j,k,1) + dimg1(j,k,1) + dimg2(j,k,1);
                    imbuf(j,k,2) = imbuf(j,k,2) + dimg1(j,k,2) + dimg2(j,k,2);
                    imbuf(j,k,3) = imbuf(j,k,3) + dimg1(j,k,3) + dimg2(j,k,3);
                
                    countbuf(j,k) = countbuf(j,k) + 2;
                end
            end
        end                
    end
    
    i = i + stepsize;
end

s = sum(sum(bwall));   


if s == 0
    bprocess = 0;
end





function [meanpixels, stdpixels, meanobjs, stdobjs] = GetMeanBWPixels(v, startframe, endframe, stepsize, cc, nstremove, nstfill, nbwarea, dthresh, xlen, ylen)

i = startframe;
endframe = endframe-stepsize;

nmax = 1000;
bwcounts = zeros(nmax, 1);
nobjs = zeros(nmax, 1);
n = 0;

while i <= endframe 
    i1 = i;
    i2 = i+stepsize;
    
    img1 = read(v,i1);
    img2 = read(v,i2);

    [bw1, idif1] = ProcessImage(img1, double(img2), double(rgb2hsv(img2)), cc, nstremove, nstfill, nbwarea, dthresh);
    [bw2, idif2] = ProcessImage(img2, double(img1),double(rgb2hsv(img1)), cc, nstremove, nstfill, nbwarea, dthresh);
     
    bw = bw1+bw2;   
    
    CC = bwconncomp(bw);
    L = labelmatrix(CC);    
    nobj = max(max(L));
    
    bwsum = sum(sum(bw));
             
    n = n + 1;
    
    if n <= nmax
        bwcounts(n) = bwsum;
        nobjs(n) = n;
    end
    
    i = i + stepsize;
end


meanpixels = mean(bwcounts(1:n));
stdpixels = std(bwcounts(1:n));

meanobjs = mean(nobjs(1:n));
stdobjs = std(nobjs(1:n));
